*Load up data file
use "`base'\Data Files and Code That Produced Them\P09 Inpatient Analysis File", replace

*These coefficients come from "P11 Estimate population for weights.sas" and are the slope of a regression of population on age in years. So convert age in months in this file
*into years then multiply by slope to get weights. Add one so it is effectively proportion and deals with 0 issue. Stata will do the rescaling. 
gen pop_all_np =  1 - 0.00926*(months_21/12)
gen pop_m      =  1 - 0.01296*(months_21/12)
gen pop_f_np   =  1 - 0.00537*(months_21/12)
 
local sufs = "all all_np f f_np m"  
local vars  = "visit  illness injury_or_alc  alcohol_any inj_by_self inj_by_oth inj_accident   "

*Create all fitted lines on either side of RD
foreach suf of local sufs {
   foreach var of local vars {  
     *Fit the left side polynomial;
      reg `var'_`suf'_r age_c age_c_sq  if age_c >= -2 & age_c < 0 [aweight=pop_`suf'], robust
      predict `var'_`suf'_left        if age_c >= -2 & age_c <= 0
     *Fit the right side polynomial;
      reg `var'_`suf'_r age_c age_c_sq  if age_c > 0 & age_c < 2 [aweight=pop_`suf'], robust
      predict `var'_`suf'_right       if age_c >= 0 & age_c < 2
   }
}
  
  
keep if age_months >= 21 - 2  & age_months < 21 + 2 
 
*Total split into two groups
gen dif = visit_all_np_r - (injury_or_alc_all_np_r + illness_all_np_r)
sum  dif
drop dif

#delimit ;
graph twoway  
 	        (scatter injury_or_alc_all_np_r age_months, mcolor(blue)  msymbol(O)  msize(.6) yaxis(1) 
			   yscale(range(50,150) axis(1)) 
			   ylabel(#3, axis(1) nogrid)) 
			(line injury_or_alc_all_np_left age_months, lwidth(thin) lcolor(blue) yaxis(1)) 
            (line injury_or_alc_all_np_right age_months, lwidth(thin) lcolor(blue) yaxis(1)) 
			  
	        (scatter illness_all_np_r age_months, mcolor(black)  msymbol(Sh)  msize(.8) yaxis(2)  
			   yscale(range(200,300) axis(2)) 
			   ylabel(#3, axis(2))) 
            (line illness_all_np_left age_months, lwidth(thin) lcolor(black) yaxis(2)) 
            (line illness_all_np_right age_months, lwidth(thin) lcolor(black) yaxis(2)) 
             , 
			 title("Figure 3: Inpatient Admission Rates", size(medlarge)) 
			 xtitle("Age at Time of Inpatient Admission", size(medlarge))  
			 ytitle("Not Injury or Intoxication", axis(2) size(medlarge)) 
			 ytitle("Injury or Intoxication", axis(1) size(medlarge)) 
			 legend(off) 
			 text(85 22 "Not Injury or Intoxication", color(black) size(medlarge)) 
			 text(105 20 "Injury or Alcohol Intoxication",color(blue) size(medlarge)) 
			 graphregion(style(none) color(gs16))  
      ;
#delimit cr
grexportpdf "`base'\Code for Figures and Tables in Paper\Figure 3.pdf"

